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Abstract. 

In the Hamiltonian approach an electron motion in a coherent packet of 
the whistler mode waves propagating along the direction of an ambient mag- 
netic field is studied. The physical processes by which these particles are ac- 
celerated to high energy are established. Equations governing a particle mo- 
tion were transformed in to a closed pair of nonlinear difference equations. 
The solutions of these equations have shown there exists the energetic thresh- 
old below that the electron motion is regular, and when the initial energy 
is above the threshold an electron moves stochastically. Particle energy spec- 
tra and pitch angle electron scattering are described by the Fokker-Planck- 
Kolmogorov equations. Calculating the stochastic diffusion of electrons due 
to a spectrum of whistler modes is presented. The parametric dependence 
of the diffusion coefficients on the plasma particle density, magnitude of wave 
field, and the strength of magnetic field is studies. It is shown that signif- 
icant pitch angle diffusion occurs for the Earth radiation belt electrons with 
energies from a few keV up to a few MeV. 
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1. Introduction 

It is suggested that whistler mode waves are responsible for electron acceleration caused 
by the interaction of these waves and Earth’s radiation belt electrons. As a result, elec- 
trons may be accelerated up to relativistic energies ( Baker et al, [1986]). Thus Horne et 
at, [1998] have identified potential whistler wave modes that are capable of resonating 
with electrons over the important energy range from 100 keV to a few MeV in differ- 
ent regions of Earth’s magnetosphere. The basic concept of energy diffusion of relativistic 
electrons resulting from resonant interaction with whistlers in the magnetosphere has been 
discussed, for example, by Walker [1993], Reeves [1996] and Summers et al., [1998]. Using 
the Chirikov overlap criterion ( Chirikov , [1979]), Karimabadi et al, [1990] considered the 
conditions at which all resonance state break up. They showed that the stochasticity oc- 
curs when the amplitude of the obliquely propagating wave is large enough. Electron pre- 
cipitation caused by chaotic motion due to coupling the bounce motion of non-relativistic 
electrons with a large-amplitude whistler wave has been studied by Faith et al, [1997]. 

The goal of work is to describe the high-energy electron motion in a coherent packet 
of whistler modes. Stochastic dynamics of charged particles in the field of a wave packet 
is one of the fundamental problem in the theory of plasma physics ( Lichtenberg et al, 
[1983]; Zaslavsky et al, [1991]). Chaotic dynamics of relativistic particles in the spectrum 
of waves is of particular interest. Thus the stochastic dynamics of relativistic electrons in 
the time-like wave packet has been discussed by Chernikov et al, [1989]. Krotov et al, 
[1998] have been shown that the stochastic heating of relativistic particles by the Langmuir 
waves in space plasmas can be regarded as a possible mechanism for the formation of the 
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energy spectrum of cosmic rays. They also studed the evolution of distribution function 
caused by the stochasticity. Nagornykh et al, [2002] developed the relativistic theory for 
the stochstic motion of electrons in the presence of obliquely propagating electrostatic 
wave. In that case, the ambient magnetic field plays an important role in randomizing 
the phase of particle with respect to the wave phase. Stochastic motion of relativistic 
electons in the whistler wave packet with application of the results to electron heating in 
the Jovian magnetoshere was studed by Khazanov et al., [2007]. 

The purpose of the present paper is to investigate the dynamics of energetic electrons 
that are confined by the Earth magnetic field and undergo bounds motion about the equa- 
torial plane. We propose mapping equation which describe the particle motion transiting 
through a wave packet in a periodic system iteratively. These equation exhibit chaotic 
motion when the particle energy exceeds a certain threshold value at given magnitude of 
the wave field. 

The paper is organized as follows. In section two we derive the canonical equation of 
motion in terms of the action-angle variables, and discuss physical processes underlying 
the results. In particular, it is shown that nonlinear nature of resonance wave-particle 
interaction leads to essentially different representation of the whistler wave field (the so- 
called space-like and time- like wave packet representations). In section three, it is shown 
that the stochastic motion of relativistic electrons can be described by a closed set of 
nonlinear difference equations. The solutions of these equations are obtained both ana- 
lytically and numerically. In section four we derive equations that describe the dynamics 
of resonant non-relativistic electrons. Solutions of these equations reveal the existence 
of an energetic threshold below which the particle can’t gain net energy from the wave 
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field. When the initial velocity is above the threshold, the electron moves stochastically 
and eventually gains a net energy. The Fokker-Planck-Kolmogorov equation is derived 
in section five, and the effects of radial drift and pitch angles scattering associated with 
particle stochastic heating are studied. Application of our results to high-energy electrons 
observed in Earth’s radiation belts is described in section six. In section seven we give 
the conclusions of our studies. 

2. Basic equations 

Let us consider a relativistic particle of charge ]e| and mass m in the wave packet of 
extraordinary electromagnetic waves propagating along an external uniform magnetic field 
of strength B. The Hamiltonian corresponding to the problem is 

H (r, p; t) = ^m 2 + (p + A) 2 , (1) 

and the canonical equations of motion are 

P = [P ,H], r — [r, H], (2) 

where p is the particle momentum, r is the position vector, A = A w + A ext the vector 
potential, the superscripts w and ext denote both the wave and external fields, and [ , ] 
stand for the Poisson brackets. 

We have employed here and throughout this paper, the frame of reference in which the 
speed of light c = 1 and charge |e| = 1. 

We denote by R the set of all real numbers. Then p e R 3 , r e R 3 , and the smooth 
manifold M = R 6 will be a canonical space of this dynamic system, and R 6 = R 3 x R 3 is 
a direct product space. 


DRAFT 


April 19, 2007, 9:36am 


DRAFT 



X - 6 


KHAZANOV ET AL.: STOCHASTIC HEATING OF RELATIVISTIC ELECTRONS 


In order to write down equations of the particle motion one must specify a coordinate 
system. We have chosen a Cartesian spatial coordinates system whose 2 axis is directed 
along the external magnetic field, the plane perpendicular to this direction is spanned by 
the orthogonal coordinates x and y. 

Making use of the connecting relations, 


B = rotA, E = — dAfdt, 
we have in the coordinate representation 


( 3 ) 


r = (x,y,z), B ex< = (0, 0, B), 

A “ = cos 0 j , 

\ k k / 

cp — zk — tuik, 

A ext = {-By,Bx, 0)/2. 


( 4 ) 

(5) 

(6) 


Here the expression for A ext is written in the axial gauge, A * is the amplitude of mode in 
the wave packet, k is the wave number, and Wjt is the dispersion equation. 

The dispersion relation for electron branch of the whistler mode waves in the cold 
magnetoplasma is written as 


k 2 /u) 2 = 1 + Up/[u(uj B - w)], 


( 7 ) 


where uj b and uj p are the gyrofrequency and electron plasma frequency, respectively. 

This equation in a long- wavelength approximation (co B uj/ujp) <C 1 reduces to 

Vp h = - uj)/u 2 , v 2 vh < 1. (8) 

We now take into account axial symmetry of the non-perturbative problem, and intro- 
duce the new variables, an action(7), and an angle (9), by a canonical transformation 
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{x,Px]y,Py) -> (0,I) ■ 

x — rcos9 , pa; = — (mru! B /2) sin 9, 

y — r sin 9, p y = (mru B /2) cos 9; (9) 

r — y2rmv^/mu B , wg — B/m, (10) 

where r is the gyroradius. 

The Hamiltonian (1) in this representation becomes 

H(z,p g ‘, 9, I; t) = H 0 (p , I) + yj 2muj R JH Q A 

■'Y^A k cos(zk + 9 - toj k ), (11) 

k 

H 0 (p,I ) = yj m 2 + p 2 + 2rnui B I. (12) 

Here we have assumed that the ratio p—.Ajm <C 1 is the small parameter of the problem, 
and retain in (11) only the leading terms. 

Associated with (11) the equations of motion are 


p=[p,H] = yj2muj B IH 0 1 • 53 kA k sin ip, (13) 

k 

i = [I, H] = yj2mu) B I H q 1 ■ ^ sin ip, (14) 

k 

z = [z, H] = pHo 1 , 9 = [9, H] = ojbttiHq 1 . (15) 

In (15) we omit the terms of the order of p 2 and introduce the definition for the phase 

ip zk + 9 — u k t. (16) 


Now the structure of wave packet, A w (t, z) = YP k A k exp [i(zk+9—tuj k )], is to be specialized. 
The most frequently used representations of wave packet are the so-called time- and space- 
like representations Zaslavsky et al., [1991], Zaslavsky, [1998]. Such wave packet may be 

exited in a plasma due to the intrinsic instabilities. In this particular physical situation, 
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we make following simplifying assumptions regarding the structure of the packet, 

k — ko + nAk, ui — uio + nAui, Ak/ko <C 1, Aw/uq < 1, n £ Z, (17) 

where k 0l ui o are the characteristic wave number and frequency, Au(Ak) is the group 
dispersion in wave spectrum, so that 

Aui = v gr Ak, v gr = dui/dk , (18) 

v gr is the group velocity, 

A k = 2 tt/L, A u) = 2-k/T, (19) 

L, T are the length- and the time- scales of the problem, Z denotes the set of all integers. 
Then we suppose that the characteristic spectral amplitude A 0 is a slowly varying function 
on t and a such that 

A 0 /uA 0 ~ 1 /tuT (< 1), VA 0 /kA 0 ~ 1/kL (<C 1), (20) 

and write down the wave packet in the form 

A w (t, z) = Aq exp(iipo) ^2 Ah exp[i(nAkz + nAojt)], (21) 

nez 

'<po ~ k 0 z + {ubtuIIq 1 - uj 0 )t. (22) 

Define the parameter, r) namely, the ratio of the particle velocity along an ambient mag- 
netic field to the group velocity, 

V = v z /v gr . (23) 

In the limit, Ak — >■ 0, rj — > 0, expression (21) can be transformed into 

A w (t , z) = A 0 exp(iV'o) S ( t / T ~ n ) ( 24 ) 

nez 
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There is the time-like representation (TLR) of a wide wave packet. The Poisson sum 
formula 

y: exp(inAcot} — y 5(f/T — n) (25) 

n£Z n£.Z 

has employed in 24. Here £(•) and the Dirac 5 function. In another limit Au: -» 0, 77 -A oo, 
we can easily show that the wave field takes the form of the space-like (SL) wave packet 

A w (t, z ) = A 0 exp(z'0 o ) y 5(z/L - n). (26) 

n€Z 

This packet represents a periodic sequence of impulses with characteristic spatial period 
L = 2ir/Ak. 

Note that both representation are often used. Thus it is established Zaslavsky et al, [1991], 
Zaslavsky, [1998]. That the TLR is available for the problem if the condition if <C 1 holds. 
The TLR of the electric field of the electrostatic waves was used by Chernikov et al., [1989] 
to derive the relativistic generalization of the standard map. On the other hand, the SLR 
has been utilized in Klimov et al, [1995] to describe the stochastic motion of relativistic 
particle in the electrostatic field of Langmuir waves, whose group velocity is small as 
known. 

3. Particle dynamics in a space-like wave packet 

We specify first the wave spectrum of a packet. Let us assume that the wave packet 
is given by (26). Thus, we will consider a relativistic electron motion in the space-like 
packet (SLP) of the whistler mode waves. In this approach dropping the subscript ”0” 
we write down the equations of motion (14) and (15) in the form: 

p — kA\J 2 ttiujbIHq 1 sin ip E - n), (27) 

n&Z 
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I = 1 Asm ip S(( — n), (28) 

neZ 

z = pHq 1 , 0 — cob'itt-Hq 1 ] (29) 

H 0 (p,I) = \Jm 2 +p 2 + 2 mu D I, (30) 

ip — uj(p, I) = kpH^ 1 + ujettiIIq 1 — to. (31) 


Here A, u and k are the magnitude, frequency, and wave number of the fundamental 
(characteristic) mode, £ = {z/L )( mod 1), L is the characteristic spacescale, 5 n = 5(£ — 
n ), <5(-) is the Dirac delta function, and Z denotes the set of all integers. 

The evident symmetries of rotation, ip — > ip ez.p(iip 0 ),ip 0 is an arbitrary constant, and 
a translation, {p, /}(£) -> {p, /}(£ + 1), allow us to represent dynamics of the system 
as certain iterative process by identification of the planes n£ and (n + 1)£. The group 
symmetry is realized first as the invariant of motion, 


p — kl — inv. 

Introducing the new variables, 


(32) 


£ z —p/m, s t — \j2mu) B I /m, e — E/m , (33) 

where E = \/m 2 + p 2 + 2mcug/ is the particle energy, and the parameter, a, 

a = 2v ph u B /u, (34) 

we write down expression (32) as 

e z - « _1 £i =0, a = ( 2u) B /u))v ph , (35) 
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where the constant of integration has chosen equals zero. In view of (35) the number of 
dimensions reduces to two. Let the variables 2 and p will be the represented pair. 

Now the explicit form of iteration system is to be found. Denote by dp, d(, dt the co- 
ordinate basis of 1-form on an extended phase space. Then, using the further result 
LdC, — pH^dt, we represent equations (27, 31) as 

de z — £~ l £ t N{A/m) sin?/> ^ 5 n d( — 0, (36) 

zEZ 

d'tp — (to — wse” 1 — k£ z £~ 1 )L(e/£ z )d( = 0, (37) 

where 8 n = 5(£ — n). Making use of the invariant of motion (35), and integrating one by 
one resulting equations, we obtain the closed set of nonlinear difference equations 

u n+ 1 = u n + 3/2 ■ Q sin w n , 
ip n+ i = ip n + N (l + 1/2 • au~l{ 3 

~ u n+l v ph\l^ + 4+1 + a«n+l ) ( m °d 2 tt), (38) 

expressed in terms of the variables u,ip, where u n+ 1 and u n are respectively the values 
of the normed momentum at times ( n + 1)C and n£, and ip n + 1 — ^ n is the phase shift 
acquired by the particle. Another deriving the above equations was given by Khazanov 
et al. , [2007]. We have employed here the following notations: 

u = £^J 2 , Q = s/a- Nw, 

N — [ kL\ , w = A/m. (39) 

Here N is the characteristic number of modes in the SL packet. The quantity w has a 
clear physical meaning: w is the dimensionless representation of the ratio of the work of 
the wave field at one wavelength, to the particle rest energy. Regarding the relationship 
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(3) between fields A w and B w , it is convenient to represent an expression for w in the 

form 


w = ab/ 2, b — B w /B. (40) 

In the relativistic limit, when the inequality e 2 > 1 is valid, the set of equations (38) goes 
over into the map: 

u n+ 1 = u n + 3/2 • Q sin ip n , 

g n ■■ (41) 

^n + 1 = ^n + {n 5/3 /Q)u~l{ 3 (mod 2 a), 

written in the perception 

u -> 7ru/ub, Q — > -kQJui „ ui, = (a 3 / 2 N 2 w /2) 3/5 . (42) 


The map g n , or rather the family of maps depending upon the parameter, becomes suitable 
for sequential analysis. So hereafter we deal with high-frequency heating of relativistic 
particles, expressed in terms of the discrete group g n , that acts on the smooth manifold. 

Denote by 

J = (43) 

o{u n ,ip n ) 

the Jacobi matrix of map (41). It is important to note that the Jacobian of (43) is equal 
to one, therefore, g n has a structure of the differentiable area-preserving map. It stands 
to reason that u and ip are the canonical pair of variables. 

We now study the behavior of this dynamic system by computational analysis. We have 
numerically integrated equations (38) and (41) for several different values of Q from 10 -3 
to 0.1. Figure 1 shows some our results computed for these equations after (10 6 — 10 4 ) 
iterations for one trajectory in the (w, ^)-phase space. The initial conditions were chosen 
in a random fashion and correspond to the region of small values of (u, ip). Figure l hen 
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Next we describe the structure of the point set. Let us consider a pair (M, g n ), where 
M is a smooth manifold and g n is the differentiable area-preserving map (41). From (41) 
follows that a pair (M, g n ) is invariable under the inversion of a point with respect to a 
circle ^(mod27r) € S, and the reflection of a point on a Oip - axis. On that basis, we 
define the following equivalence rules: 


(-7T, u) ~ (7T, «), (</>, -7T) - (V», 7T), 


(44) 


where ~ stands for the equivalence sign. 

The identification of these points is indicated by the arrows in Figure 1. 

The equivalence rules (44) allow us to represent phase space of the system as the space 
having the topology of a torus, T 2 = S x S, ^(mod27r) e S,u 6 S(Kosniowski, [1980]). 
Now we discuss local topology of the manifold considering a Jacobi matrix given by (43). 
Denote through Ai and A2 the eigenvalues of the matrix J. Recalling that local structure 
is determined by neighbors of fixed points, we have from (43) and (39) 


det J = Ai- A2 — 1, 


(45) 


tr J = Ai + A 2 = 2 + (tt 5 u 5 ) 1//3 , (46) 

where det J and tr J denote the determinant and the trace of this matrix, respectively. 
As pointed out in Arnold et al. , [1968] the condition 

| tr J | = 3 (47) 


corresponds to a topological modification of a phase space, and its validity implies that 
the manifold have topology of a hyperbolic torus. Thereupon from (47) subject to (46) 
we find 


Ai — 


3 + x/5 


A 9 


3 — yjh 


2 2 

April 19, 2007, 9:36am 
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Since Ai and A 2 are the dynamic characteristics such that A 2 > 0, Ai > 1, and the ratio 
A]./ A 2 is an irrational number, therefore the map (ip n ,u n ) — g n (ip 0 ,w 0 ), (ipo,u Q ) is an 

initial point of phase curve, forms a stochastic phase flow with the mean rate of a loss of 

infomation, K, 

K — In Ai, (49) 

K is the Kolmogorov entropy. 

Present conditions (46) and (47) ensure that the relation 

\ub\ — tv (50) 

take place, and it determines the upper bound of {u}. 

Like that the set {u,ip} is a compact, whose structural stability is determined by the 
fractal dimension, df 

df = 1 — In Ai/ln A 2 , df = 2. (51) 

We call any compact a probabilistic fractal if its topological dimension is less than df and 
K > 0 Mandelbrot, [1982]. In our case both these conditions are valid. 

As seen in Figure 1 for a given values of Q the stochastic region extends to values of u 
predicted by above equation (50). We know, that u depends on the parameters of this 
problem as given by (39). By that condition (50) determines in itself an equivalence class 
in the Q-parametric space. 

Considering (39,42) and (35), from (50) the following equation results 

/I \ 2 / 5 

e t = a^6JV 2 j , (52) 
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which determines the upper value of the energy spectrum for the SLR case. Note the 
dependence £b on the driving field b is weak enough. This agrees to the numerical solution. 

Finally, it should be noted that the phase flow, (M, g n ) is structural stable and typical 
because df and K be an invariants, and sup{c} smoothly depends on the driving amplitude 
6, which is the controlling parameter of the system. Phase flows with such properties are 
said to be the stochastic (strange) attractors. Certainly the flow (M, g n ) is such attractor 
as was to be proved. 

At last, the achieved fractal measure df — 2 infers that points of phase curve evenly fill all 
obtainable phase space, or in other words, all states of our dynamic system are equivalent. 

This issue is supported by simulation in Figure 2 on which is represented a joint probability 

density, p(ip,u). The following algorithm was used for computing p. We partition all the Figure 2 hen 

phase space on the identical cells with the mesh size AipAu — (27r/40)x(27r/30). Like 

that p(?/>, u) is proportional to the number of phase points in the element of phase space, 

AipAu. Figure 2 reveals that diffusion in the phase angle is very fast. Significant phase 
angle diffusion occurs on timescales order of tens of T, where T is the step of one iteration. 

On these timescales, the change in u is smaller (the numerical calculations indicate the 
characteristic time for establishing the uniform distribution in u is proportional to T/Q 2 
at Q <C 1), therefore, the variable u is a slow varying coordinate on the strange attractor. 

Statistical aspect of the problem will be studied in more detail in Sec.5. 

4. Non-relativistic electron motion in a time-like wave packet 

We now discussed the dynamics of non-relativistic electrons in the wave packet of the 
whistler mode. In the case, the condition v z /v gr < 1 is valid, therefore by means of (24) 
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the wave field of packet may be approximated as 

A w (z, t) — A 0 cos(koz — u B e~ l t — uot) ^ 6(t — n). (53) 

n 

Here r = t/T, T is the timescale of the problem. 

Dropping the subscript ”0” in (53) and substituting this expression in (10) we write down 
the Hamiltonian as 23 


H(z,p; 9, I; t) = H 0 (p,I) + ^UEl/mAcosip 

■ J2 6 ( r - n )> 

(54) 

n£Z 

H 0 = p 2 /{2m) + oo B I, mv^/2 — u B I, 

(55) 

ip = zk + 9 — tu. 

(56) 

The equations of motion associated with (54) are 

p = k^/2u B I /mA simp ^ 5(r — n ), 

(57) 

i = \J2 uj b I /mA sin^y^ 5(r — n) 

(58) 

z = v z =p/m, 9 = ui b . 

(59) 


In (59) we retain only the leading terms on condition that A/mv t is a small parameter of 
the problem. 

Since the whistler wave frequency u is smaller than the electron cyclotron frequency u B , 
electrons in general must move in the opposite direction as waves in order to match the 
resonance condition 

kv z + ujb — w — 0. (60) 
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In view of this fact and also taking into account the relation ujgl = mv 2 /2, we write the 
invariant of Hamiltonian flow (58, 59 ) as 

vf - a(v r ~ K|), a - 2(u B /u)v ph (61) 

where v r = v P h{u B — oj)/lo is the resonance speed of an electron. 

To understand the physical picture of the stochastic motion we use the overlap criterion 
Chirikov, [1979]. In the broad wave spectrum case given by (25), the condition (60) 
describes the family of the resonance states. First varying this equation we evaluate the 
interval between an adjacent resonance states in the velocity space as kAv z ~ Aoj = T~ l , 
and 

Av z ~ v ph /N, N = [wT], (62) 

where N is the characteristic number of modes in the TL-packet. 

Then it avails oneself of equation (57) to estimate the width, Sv z , of resonance state in 
the velocity space, 

6v ~ abNv t /2v ph , (63) 

where we have again used the notation w — A/m — ab/2, b — B w /B. 

Now from the overlap criterion, 5v > Av , and the two above expressions the following 
relation results 

v t >v c = v ph oj/uj B N 2 b. (64) 

Note the term in the RHS of equation (57) corresponds to the Lorentz force caused by 
spiralling an electron in the magnetic field of wave. Consequently we conclude from above 
the equation that the stochasticity occurs when the Lorentz force acting of electron in the 
axial direction exceeds some defined value. 

DRAFT April 19, 2007, 9:36am DRAFT 



X - 18 KHAZANOV ET AL.: STOCHASTIC HEATING OF RELATIVISTIC ELECTRONS 

Now we may be able to estimate the value of wave field at which an electron motion 

becomes stochastic. We substitute the value of v c from (64) in (61) at v z — 0 to obtain 

b > b c - (v ph /2v r ) 1/2 {uj/cu B ) 3/2 N~ 2 . (65) 

Now we turn to the equations of motion. We use (61) and the expression mv\j2 — u B I 
once more to represent the equations of motion in the form 

v t = ~u B ab sin ip ^ 5(t -ra), (66) 

1 nez 

'ip — u(v t , w) = u B - uj - kv r + (k/a)v 2 . (67) 

By integrating these equations can be transformed into a map, G n , 

s~m . ^n+1 ”t" QsilCl'lp n: -(w G JJ C R\ 0<«< 1}, IrsI 

^n +1 = 'ipn + sul +1 (mod 27 t), ^ ' 

where subscript n refers to values taken at time t = nT , and the new variable, u, has 
introduced by the relation 

u = v t /v o, v 0 = (y.Vgr . (69) 

Here Vq is the all perpendicular (to B ) speed of electron, and the new parameters are 

Q = ^Nabcos/voUJ, N = [wT], s = Nvl/av ph . (70) 

We have numerically integrated equations (68) for values of Q from 0.0005 to 0.01. Our 
results are shown in Figure 3. 

The overall picture of the phase space is quite different for u < u c and u > u c . In the 
first case, the motion is regular. The figures indicate the existence of a threshold for the 
initial particle velocity above which the trajectory becomes chaotic and the extent of the 
stochastic region increases linearly with Q. 


Figure 3 
here 
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Thus the Jacobian of the matrix of Q n , det J — 1, therefore G n is the measure-preserving 
map, and ip, u are the symplectic pair on the smooth manifold, M : U x S. Thereby we 
apply the condition (50) to (68), to find the expression 

u c = znf{u} = (2 Qs)~ l , (71) 

that well determines the lower bound inf{u} of the stochastic set. 

Substituting (69) and (70) in (71) we get 

v c = v ph uj/u) B N 2 b , (72) 

which is in good agreement with the numerical solutions and results of qualitative analysis. 

So we resume that nonlinear electron acceleration by a wave packet of the whistler mode 
waves is always a stochastic process. 

5. Diffusion evolution 

Particle dynamics in random electromagnetic fields is known to be described by the 
quasilinear theory (QLT) ( Sagdeev et al, [1969]). The QLT approach in particular, has 
been employed to deal with diffusion of electrons in the turbulent field of whistler waves 
packets ( Kennel et al, [1966]), ( Summers et al, [2004]). In coherent electromagnetic 
fields, on the other hand, the particle dynamics are not described by these theories. The 
nature of diffusion in this case is described by stochastic dynamics of particles, when 
the motion along stochastic trajectories gives rise to the so-called deterministic diffusion 
( Lichtenberg et al, [1983]). General mathematical and physical aspects of this problem 
have been discussed, for instance, in reviews ( Arnold et al, [1968], Zaslavsky et al, [1991]). 

The purpose of the present paper is to investigate electron dynamics in a coherent 
packets of whistler mode waves. We have employed here the method in which dynamics 
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of the phase variables on stochastic attractor is included in a Fokker-Planck-Kolmogorov 

(FPK) equation ( Zaslavsky et al., [1991]). 

In the present section, we shall exploit the canonical Hamiltonian structure that has 
been developed above. Thus, the drift kinetic equation follows at once since the canonical 
structure equations of motion (28, 29) have been established. The distribution function 
(probability density) / (w; t) obeys the Fokker-Planck-Kolmogorov (FPK) equation 


df{u-,t ) = ld_ D dJ_ 
dt 2 du du' 


which holds if Q <C 1 ( Zaslavsky et al, [1991]). 

Here D is the conventional diffusion coefficient in phase space, 


(73) 


D =< (u„ +1 - u„) 2 > T *, 


(74) 


in which (w n+ i — u n ) is substituted from (41), and <• > denotes the phase average, T is the 
timescale of mapping (41). The function f(u : t ) belongs to the space of all differentiable 
functions supported in [— 7T, 7t], and the functional 


f f(u, t)du — 1 

J — 7T 


(75) 


is the condition of normalization. 

First, by means of (41) we calculate using (74) the diffusion coefficient 

D = 9Q 2 /8T. 


(76) 


Then, making use of the result (50) proved above along with (73) and (76), we evaluate 
the characteristic time for redistribution u over the spectrum 


t d ~ ul/D = T(4 tt/3Q) 2 . 


(77) 


DRAFT 


April 19, 2007, 9:36am 


DRAFT 



KHAZANOV ET AL.: STOCHASTIC HEATING OF RELATIVISTIC ELECTRONS X - 21 
A objective of this study is to determine the time-independent distribution function and 
the rate of heating. This requires a solution of the FPK equation together with a nor- 
malization (75). Thus a solution of the FPK equation along with the boundary condition 
/(— 7r) = /( 7r) in the limiting case t > may be given in the form of the uniform 
distribution 

f(u) = (27T)- 1 . (78) 

Now, it becomes relevant to determine how the system evolves in time at t < We 
exploit the FPK equation with f(u ) and its derivative df /du vanishing at the boundary. 
We introduce the moment < u 2 >= J* v duu 2 f (u) , multiply equation (73) by w 2 , and 
integrate the resulting equation over u to obtain 

d <u 2 > /dt — D. (79) 

We know, that the variables u and £ are related by (42). This allows us to attach all 
possible states of e a probabilistic measure, namely, the probability density, which 

is associated with f(u,t ) via the measure-preserving point transformation 

/(M) - f(u,t)(du/d(s)). (80) 

Representing our results we start with the SLR case. In the case the variables u and £ 
are associated by relation (42). Then equation (39) along with the normalization 

/ f(e) de = 1 

J 0 

allow us to derive the time-independent distribution 

f(e)de = e^de, {s G R + |1 < e < e b }, (81) 
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and the FPK equation describing the evolution of an energetic spectrum in the allowed 

range of particle energies 


d 

dt 


f(e,t ) = 
D(e) = 


a 3 N 2 b 2 
IQsT ' 


(82) 

(83) 


One may observe that there is an explicit dependence of D on e~ l . Consequently the 
heating rate of electron decreases with its energy as e“ 2 , namely, 


a 3 N 2 b 2 
487V ’ 


1 < £ < £b- 


(84) 


The restriction on e in (84) is needed because the energetic spectrum is bounded above 


by £b- 

Then from (82), (83) and (73) we find the characteristic time to establish the distribution 
(81) 

t d = 16T(Af/646 2 ) 2/5 . (85) 


Now via (81) we evaluate the means 


< e >= (l/2)e s , < e 1 >= (1/3 )s 2 b , 


( 86 ) 


and the relative standard deviation, 

V< e 2 > - < e > 2 

= 0.5. 

> 

The distribution function f(u, t) in the non-relativistic case also obeys the FPK equation 
as a consequence of the fact that the equations of motion, expressed in coordinates ( u , ip) 
reflect the underlying canonical Hamiltonian structure of the equations (61). By reason 
of (71), we write the normalization as 




April- 


JjH 


du — 1. 

007, 9:36am 
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Then by means of (70) along with (71) and (70), we calculate by formulas (74) and (77) 
the diffusion coefficient 

D = Q 2 /2T = a 2 N 2 b 2 u 2 B /8Tvluj 2 , (89) 

and the characteristic time for redistribution u over the spectrum 

t d = AT /Q 2 = lQTv 2 0 co 2 /N 2 a 2 b 2 u 2 B . (90) 

Next using the normalization for f(u,t ) in the FPK equation, it follows when t > t d the 
distribution function can be given by 

f(u) du, = (1 — u c )~ l du. (91) 

This means that the random variable u is evenly distributed on [u c , 1] . 

Now from the transformation f(u) du = f(v) dv and the relation u — v t /v d we obtain 

f(vt) = v d l . (92) 

Taking account of (88), by (92) we derive the quadratic mean of v t 

< v t> = 7 f(v t )v 2 dv t = v 2 J 3 

J o 



= 4/3-v% h (uB-u))/w, 

(93) 

the mean value of \v z \ 

and the mean velocity 



< \v z \ >= 2/3 ■ v ph (u B - u)f(jJ B , 

(94) 


< v z >= 0. 

(95) 

Then we define by 




T t =< mv 2 / 2 >, T z —< mv 2 /2 >, 

(96) 
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to find the ratio 


T t /T z — 5a/8 v gr — 5/8 • m|/a j(uj b — u), 


( 97 ) 


describing the anisotropy of distribution (61) in the (v z ,v t ) phase space. Anisotropic 
distribution function (61) with its derivative 


— = -\\l a l( v sr- I v z I) sgnv z , 


(98) 


dv z 2 

that tends to infinity as \v z \ — > v gr , describe the so-called pancake distribution in the 
(v z , v t ) phase space. 

Again it is useful to evaluate the energy distribution function, f(e), for the domain 
of non-relativistic energies. Considering (61), we note that the particle energy E is a 
specified function of v t . Thus the measure-preserving transformation f{e)de = f(v t )dv t 
determines this problem completely, subject to appropriate boundary conditions. To that 
end, we have to find the lower (e c ) and upper (so) boundaries of an energetic spectrum. 
At first, from (61) we evaluate threshold value of particle speed 


V th | hz |c Vgr V c /Oi. 


(99) 


Now it is clear that s c and Sq are given by 


= 1/2 ■ (ty - K/o f + 


e 0 = Wj/2 = aVg r /2. 


(100) 

( 101 ) 


The threshold energy e c is a function of the magnitude of wave field, b, and, in view of 
(71) s c hinges on b as follows: 


e c — ^ oo as b 4, as b — ■ y 0, 


{s c - V 2 /2) -4- 0 as b 


-2 


( 102 ) 

(103) 
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The corresponding solution for f(e ) provided that 

re o 

/ f(e) de=l 

J e c 

can be expressed in the form 


1(e) de- 1/2 ■ [(e 0 - e c )(e - e c )] 1/2 de, e a > e> e c . 


(104) 


This solution relates directly to the behavior of the system near the order-chaos bifurcation 
transition. The function f(e) undergoes a sudden change at e = e c , and in the region 
of regular motion (e - < e c ,e = const) has the form of the Dirac (5— function, f(e ) = 
d(e — const). 

Thus f(e) describes the density of states in an energy space, the mean particle energy 
and relative standard deviation can be calculated by 


f£ 0 

/ ef{s) de = 1/3(f 0 - e c ) + £ 

J S c 

\/< £ 2 > - <£ > 2 




= 0.4. 


£ z > 


(105) 

(106) 


The latter indicates the high level of fluctuations in the energetic spectrum of an electron. 


Evolution of a means at t < t d is governed by (79), where D is now given by (90). 


Equation (79) in the explicit form may be written as 

~7 < vl >— a 2 b 2 N 2 u! 2 B /8Tu 2 . (107) 

Cl L 

Then we have respect to the invariant of motion along with (107), and define by /j, a = 
v 2 /2B the first adiabatic invariant to find the rate of a change of fj, a 

fi a — a 2 b 2 N 2 u%/lQTBui 2 , (108) 
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and the heating rate 

e = D e ,D s = a 2 b 2 N 2 u 2 B /16Tco 2 . (109) 

This result is nontrivial, because chaotic motion in the (u, ip) phase space leads to impor- 
tant and easily observable macroscopic effects such as the stochastic heating of plasma 
particles. 

Now we describe the effects associated with stochastic heating of high-energy particles. 
First, we discuss the pitch angle distribution in a range of non-relativistic energies. 
Denote by Xp the pitch angle, and let x he its complementary angle. Then we define 

tan x — v z /v t , x € (-7t/2, tt/2), (110) 

and use equations (61) to represent tanx as a function of u(u = v t /v b ), namely, 

tanx — (vgr/a) 1 ^ 2 • (1 — u 2 )/u. (Ill) 

Then via (91) and the measure-preserving point transformation f(u)du = f{x)dx , we 
find the angle distribution 

fix) = (a/v) 1 ' 2 ("7(1 +“ 2 ) 

+ (v gr /a)(l — u 2 ) 2 /(l + u 2 )) . (112) 

Equation (112) determines the angle distribution function in a parametric form. Nonlinear 
transformation (111) allows us to express (112) as a function of a proper variable X- Note 
that /(x) is invariable under transformation 

u y w, x -x; fix) = f(-x), (ns) 
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therefore fix) is a symmetric function on ( — tt/2, 7t/ 2). However, the explicit expression 
for f(x) has a very complicated form. The following asymptotic formulas are valid: 

fix) ~ 1/2 • (a/v) 1/2 (l - («/v) V2tan x) 

as x 0, (114) 

fix) ~ (^r-/«) 1/2 ( 1 + tan -1 x), asx~>±7r/2, (115) 

fix) ~ iv gr /a) 1/2 il + 2ie/s c -l)) 

asefe c , x ±7t/2, (116) 

/(x) w 1/2 • (v/ a ) _1/2 ~ i) 

as e f eo, x 0- (117) 

where e c and e 0 are given by (100) and (101). 

Making use of these equations, we conclude that the function fix) governing the dis- 
tribution of angles is a convex symmetric function of x? having the maximum at 
X = 0 ix P = 7t/ 2). The result indicates also that particles are equally likely to be 
scattered in the direction with or against the wave. To this may be added that the 
given function is an increasing function of the particle energy on (e c , e 0 ), and it hinges 
asymptotically on £ as (e — e c ) near the energetic threshold, and as Jef£c — 1 at e > 2e c . 

Let us discuss particle scattering in a range of relativistic energies. Because ^4/m <C 1, 
we write 

tanx = v z /v t ttPz/Pt ~ e s /e u X € (“V 2 . k/2). (118) 

at energies above 1 MeV. 

Taking account of the invariant of motion (35) we derive the following dependence 


(119) 
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Using (81) along with the measure- preserving transformation f(x) dx = f{e) de yields the 

time-independent pitch angle distribution /(x), namely, 

f(x)dx = 2tan” 2 y m -cos _2 x- |tany|dx, 

X e {-Xm, Xm), tanxm = {£b/a) lf2 i (120) 

and its dependence on the particle energy 

/(x(c)) = 2^ 1 (ae) 1/2 (l + eof 1 ). (121) 

The function given by (120) is the concave symmetric function of x, which has a maximum 
at x — ?r/2 {x p — 0), and its derivatives tend to infinity as x — > ±7t/2. This function 
describes the so-called U-like distribution. 

Since the magnetic field is constant, the gyroradius is a direct measure of the perpen- 
dicular electron velocity. Thus the stochastic heating will be accompanied by a radial 
drift of particles in space. Indeed, in view of the relation r = v t /u>B and (107), we get the 
following expression 


d < r 2 > 9 d < v? > _ 

= uZ 2 — = D u 

dt dt 

(122) 

D t = a 2 b 2 N 2 /1 §Tu 2 b , e > e Cl 

(123) 


where D t is the coefficient of collisionless diffusion across the ambient magnetic field. Then 
we evaluate the radial drift in a range of relativistic particle energies. 

First via (9) and (35) we get the relation . 

r = y/aiujg 1 . (124) 

Now we write down the diffusion coefficient as 

N a A N 2 b 2 o 

A(£)= 48F4' £ ’ (125) 
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which immediately follows from (84) and above result. 

Also relation (124) allows us to represent D t as a function on r, 


Di{r) 


a 4 N 2 b 2 fr 0 


(^j , r o = \fa/u) B , r b = y/aei/wBr-ro <r <r b . (126) 


48Tcu| 

Equations (123, 125) and (126) are correct if the inequality r/L <C 1 holds. In view of 
equation (124), this inequality can be written in the form r/L — y/aelc o/NusVph -C 1. 
For typical parameters this condition is satisfied trivially. 


6. Application 

Enhanced convection electric fields associated with solar wind streams provide the prin- 
ciple mechanism for the intensification of ring current (10 — 100) keV flux, and also leads 
to the excitation of whistler mode waves Hastings et al, [1996]. The zone of most intense 
wave activity is spatially localized because of the decrease in resonant energy and wave 
guiding by strong density gradients associated with the plasmapause Wolf [1995]. Typical 
whistler wave amplitudes are in the range (10 — 100) pT, but occasionally wave amplitude 
approaches 1 nT. The wave magnitudes of outer zone chorus emissions are usually small 
enough, resulting in weak diffusion scattering. Consequently, the electron heating should 
occur gradually over many drift orbits. 

To check the performance of the method, at first the parameters of the problem are 
chosen to reflect a typical whistler wave interacting with electrons at locations outside 
the plasmapause, where the electron cyclotron frequency u B = 6 • 10 4 s _1 , the electron 
plasma frequency cu p — 12 • 10 4 s -1 , and the wave frequency cu — 3 • 10 4 s _1 . Thus 
the dispersion relation (8) corresponds, undoubtedly to the problem. By (8) we find 
the phase velocity, v p h = 0.25, and the group velocity, v gr — 2 • v p h (u) B — u>/u) B ), v gr ~ 
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0.25. The theory is based on three parameters, namely, a,N,b. First we evaluate a = 

2 v p h^Bl^i a = 1-0. Typical magnitudes of whistler mode waves are in a range of 10 -6 — 
10~ 7 G, consequently the parameter b is of order 10~ 4 — 10” 3 . We set b = 1.3 • 10 -4 
below. To evaluate the parameter .N — ojT for the time- like wave packet, we need to 
estimate the transit time of the wave packet through the domain of resonant interaction 
(the zone of intense wave activity). Thus the characteristic size of this region is of the 
same order of magnitude as Earth’s radius, we have accepted L — 5 • 10 8 cm. By that 
the value of T, T — L/v gr: is about 7 • 10~ 2 s, and N = 2 • 10 3 . In the following we 
will use these values to obtain an estimate of some quantities. First we calculate the 
threshold value of electron energy, E c , at a given value of b. Applying (71) to (100) yields 
an estimate of E c , E c = mc 2 e c , E c ~ 16 keV. Then we calculate by (52) the upper value 
of the energetic spectrum, ~ 3.5 MeV. Considering (123, 109), we evaluate the rate 
of stochastic diffusion across the ambient magnetic field, D t ~ 3.5 • 10 10 cm 2 s _1 , and the 
heating rate, D E — D s mc D E — 70 keVs -1 . In the relativistic region of particle energies, 
the expressions (125, 84) yield D t ~ 1.1 ■ 10 10 (m/E) 2 cm 2 s _1 , D E — 25(m/E) 2 keVs -1 . 
Note that D E as well as D t are a decreasing function of E in a relativistic range of the 
energy spectrum, while D E and D t in a non-relativistic domain of particle energy are 
independent on E. The following dependencies, u) p oc y/n^ : v ph oc Bj^fnl and aN cc B, 
allow us to express D E as a function of wave amplitude B w , strength of an ambient 
magnetic field B, and plasma particle density n e and arrange the resulting expression as 
the scaling law D E — const. (B / y/n^)(B w ) 2 . Using above result in (123), we reveal that 
D t hinges on B as B~ x . 
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Now the effect of pitch angle scattering of relativistic electrons should be evaluated. From 
(120) at £(, — 7 and a = 1 follows that in this case whistler waves can diffuse electrons in 
a cone with the vertex angle, Xv = 2 Xb, which is about 140°. The equation (121) indicates 
a degree of scattering anisotropy, deg A = f(x)( £ = £ b)/f{x){ £ = 1), increases with e up 
to deg A ~ 20. Now we respect to the pitch angle scattering of non- relativistic electrons. 
The results (116) and (117) show that the wave packet to effectively scatter electrons in 
pitch angle leading to establishing the distribution peaked at Xp — 90° with the degree of 
anisotropy deg A — 2. 

Now to see that our approach correctly interprets the facts, it is necessary to find the 
range of validity of our results. Thus the basic equations are applicable provided kL >> 1, 
and t obT 1. Then it is necessary to establish that the characteristic time for the energy 

the distribution over the spectrum is small as compared with the period of an intense wave 
activity. According to (90) and (85) the energization time for non-relativistic electrons is 
about 7 s, while the time for establishing the energy distribution in a range of relativistic 
energies is approximately 0.5 h. Note that the efficiency is of this process is high enough, 
because observations show that the power of whistler-mode waves is enhanced at substorm 
injection and decays over a period of hours or over a period of a few days during the storm 
recovery ( Baker et al., [1986], Horne et al., [2003]). 

Now the threshold value of the wave field given by (65) may be evaluated as b c ~ 10 -7 . 
The equations of motion (38) in the region of relativistic energies are obtained in the 
faithful representation, and these include the ultra-relativistic approximation, e 2 ;§> 1, 
that has been employed in deriving (41). According to (35), this inequality may be 
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written in the form 

2(u) B /u)v P h£~ 1 < 1. (127) 

Owing to (52) the latter can be transformed into the requiring 

b = B w /B > 47T 2 (~HT 6 ). (128) 

These conditions are fairly easy to check. Recalling direct amplitude measurements 
(. Nagano et al, [1996]) and the extremes of spectral density measurements ( Parrot et 
al., [1994]) yield ( B w jB ) = (10 -4 — 10~ 3 ), we conclude that (128) trivially satisfied, 
therefore, the lower- and the upper limits for the energetic spectrum are correct. 

Note that in accordance with Faith et al., [1997] the magnitude of wave field required to 
induce chaos by a single wave is about 10 -1 , which is quite larger compared with those 
actually observed. 

It should be noted that the geometry of an ambient magnetic field does not play a role in 
the given problem, thus the lengthscale of an interaction region is typically smaller than 
the travelling path, i.e. an ambient magnetic field is always locally an uniform field. To 
this may be added that the analytic model describing the motion of a particle out of the 
interaction region must include the radial diffusion of electrons conserving the first and 
second adiabatic invariants. 

These results demonstrate the possibility of electron heating over a wide energy range 
between > 10 keV and a few MeV. This is in reasonable agreement with the experimental 
observations. 
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7. Summary 

A canonical Hamiltonian approach has been employed to deal with the whistler wave- 
electron interaction and the stochastic heating of high-energy electrons in magnetized 
plasmas. The time-like and space-like wave packet representations have used in deriving 
the equations of motion for non-relativistic and respectively relativistic electrons. It is 
shown that phase dynamics of the system is realized on stochastic attractor, all means on 
this attractor are stable and irrelevant to any initial conditions. 

The results afford a basis for conclusions. 

Solutions of these equations have revealed the existence of an energetic threshold below 
which electron motion is regular. This threshold can be expressed in terms of the normal- 
ized magnitude of wave field, b, as g: c = const.b~ A as b -> 0, and (e c — v 2 gr /2 ) = const.b ~ 2 
as b — > oo. When the initial energy is greater than the energetic threshold, an electron 
moves stochastically. Thus, nonlinear electron acceleration by a whistler wave packet is 
always a stochastic process. Chaotic motion gives rise to diffusion in energy, and leads 
to establishing an energy spectrum, bounded above by e &. The upper value of energy 
spectrum has a weak dependence on b (oc b 2 / 5 ). The energy spectrum in a range of non- 
relativistic electron energies obeys the law /(e) = const./ (e — e c ) which describes the 
order-chaos bifurcation transition at e = e c . While the relativistic electron distribution 
obeys the uniform distribution. 

The obtained energetic spectra have used for evaluating pitch angle distribution functions 
over different energies. Whistler mode waves can scatter non-relativistic electrons with 
pitch angles of up to 90° and are responsible for the formation of the pancake distribution 
in the (v z , v t . ) phase space. Relativistic electrons are mainly scattered by this waves 
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almost along the direction of an external magnetic field. The high-energy electrons with 

chaotic motion are equally likely to scatter in the direction with or against the wave, and 
there exists a certain probability that some particles will be scattered into a loss cone and 
will precipitate into the polar regions along magnetic field lines. Chaotic motion of an 
electron in a constant and homogeneous magnetic field is accompanied by radial drift. This 
effect must has associated with stochastic heating of the particle. High-energy electrons 
observed in Earth’s radiation belts can be regarded as tails of the electron distributions in 
magnetosphere. According to our analysis, for a spectrum of a whistler modes interacting 
nonlinearly with Earth’s radiation belts electrons, substantial energization can occur in a 
range of energy from a few keV up to a few MeV. Thus stochastic processes occuring in 
the plasma can play an important role in the evolution of the energy spectra of radiation 
belts electrons. 
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Figure 1 . The phase space of the map g n . (a) One single trajectory of length 10 6 for 
Q — O-OItt. The trajectory was started from the point Uq — 10 -3 , ipo = 10 -4 . The identification 
of points is indicated by the arrows. (b) Close to zero, approximation g n loses its validity. 
Therefore in the vicinity of u = 0 we must use the solution of the original equations (52) and 
(38) expressed in terms of the variables of u, ip given by (42). The figure shows the boundary of 
the chaotic region is well approximated by the condition (50). The parameters are the the same 
as in Figure la. 


Figure 2. Joint distribution (density of an imaging points) F(ip, u ) obtained numerically from 
(41) with the parameters as in Figure 1. 


Figure 3. The phase space of the map G n . (a) A single trajectory of length 10 6 for Q = 10 -3 , 
5 — 5 • 10 3 , (w 0 , ip 0 ) — (0.01, 0.0001). (b) A phase curve for the same values of Q and s as in (a). 
The trajectory was started from the point (0.11, 0.0001). (c) One trajectory of length 7 ■ 10 6 for 
Q = 0.5 • 10" 3 s = 510 -3 , (uq, ipo) = (0.21, 0.0001) 
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